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1  .  INTRODUCTION 


Flow  separation  introduces  major  complications  for  the  designer  of 
high-speed  flight  vehicles.  Separation  can  result  in  increased  drag,  loss  of 
lift,  buffeting,  reduced  control  surface  effectiveness,  and  very  high  heating 
in  the  reattachment  zone.  These  effects  are  most  pronounced  at  high  Reynolds 
numbers,  where  the  boundary  layer  is  turbulent.  Computing  a  separated  flow 
is,  of  course,  a  much  more  difficult  task  than  computing  an  attached  boundary 
layer.  Modeling  the  turbulence  behavior  is  the  other  major  problem. 

Separated  two-dimensional  (2-D)  flows  have  been  studied  extensively.  The 
more  recent  works  have  concentrated  on  numerical  solutions  to  the  Navier- 
Stokes  equations.  Shang  and  Hankey^  used  a  parabolized  Navier-Stokes  (PNS) 
method  for  shock-separated  flows;  they  found  that  a  lag  method  is  required  to 
achieve  decent  agreement  with  observed  pressure  distributions  in  separation 
zones.  In  a  series  of  papers,  Horstman^  and  co-workers  investigated  several 
turbulence  models  in  Navier-Stokes  solutions  for  shock-separated  flows.  The 
turbulence  models  include  an  algebraic  eddy  viscosity  model,  a  one-equation 
model  consisting  of  an  additional  partial  differential  equation  for  the  turbu¬ 
lent  kinetic  energy,  and  a  two-equation  model  (pde's  for  turbulent  kinetic 
energy  and  length  scale).  They  found  that  the  more  sophisticated  models 
(i.e.,  the  one-  and  two-equation  turbulence  models)  yield  somewhat  better 
agreement  with  detailed  measurements  of  shock-separated  flows.  However,  none 
of  the  models  were  entirely  satisfactory. 

This  report  summarizes  the  results  of  a  study  to  develop  an  improved 
engineering  method  for  predicting  separated  supersonic  turbulent  flows.  The 
approach  is  based  on  the  use  of  an  inverse  boundary  layer  method,  coupled  with 
a  five-equation  Reynolds  stress  turbulence  model.  The  inverse  method  promises 
to  permit  meaningful  solutions  to  the  boundary  layer  equations  for  separated 
flows,  with  much  less  numerical  effort  than  would  be  required  to  solve  the 
Navier-Stokes  equations.  The  Reynolds  stress  modeling  should  be  able  to 
account  for  the  important  effects  of  adverse  pressure  gradients.  Pressure 


gradients  significantly  alter  the  distribution  of  energy  between  turbulence 
components,  an  effect  that  is  not  readily  incorporated  into  simpler  one-  or 
two-equation  models.  The  goal,  then,  is  to  obtain  an  accurate  description  of 
turbulence  behavior  in  supersonic  separated  flows,  at  a  computational  effort 
appropriate  to  engineering  applications. 

All  of  the  discussion  here  will  be  restricted  to  2-D  separated  flows. 
Three-dimensional  separation  may  be  more  important  for  real  applications,  and 
the  2-D  idealization  may  represent  a  singular  limit  to  the  more  general  3-D 
situation.  Nevertheless,  there  is  good  experimental  data  on  2-D  separation 
and  it  seems  appropriate  to  attempt  to  understand  the  simpler  case  first.  As 
will  be  discussed,  there  are  several  unresolved  challenges  to  the  2-D  case. 
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2.  SECOND-ORDER  CLOSURE  MODEL 


The  development  of  the  model  used  in  this  study  has  been  described  in 
detail  elsewhere,2'4  and  only  the  results  shall  be  presented  here.  Our  treat¬ 
ment  of  closure  draws  upon  aspects  of  various  previous  work,  most  notably 
Rotta's  treatment  of  low  Reynolds  number  effects®  and  the  description  by 
Hanjalic  and  Launder®  of  the  triple  fluctuation  and  pressure  fluctuation 
terms.  Wherever  possible,  closure  approximations  have  been  evaluated  against 
basic  laboratory  experiments  (e.g. ,  grid  turbulence),  the  types  of  experiments 
being  chosen  to  attempt  to  isolate  individual  terms. 

The  formulation  accounts  for  both  mean  and  fluctuating  velocity  and  tem¬ 
perature  quantities.  The  dependent  velocity  variables  are  the  mean  velocity 
vector  U ,  the  Reynolds  stress  tensor  u^u  j ,  and  the  isotropic  dissipation  rate 
<f ;  under  the  boundary  layer  approximation,  this  set  of  variables  reduces  to 
U,V,uz,  vz,  wz,  uv  ,  and  <5.  In  practice,  it  is  convenient  to  replace 
u' 2 ,  v' 2 ,  w' 2  by  the  kinetic  energy  q2  ■  (u ' 2+  v‘ 2  +  w' 2)/2  and  two  measures 
of  the  degree  of  anisotropy  }  *  u* 2  -  2/3  q2 ,  S22  “  v ' 2  -  2/3  q2 . 

For  steady  flow,  the  governing  equations  include  continuity: 


pUj  -  0 


the  mean  momentum  equation: 


„  3U  3p  3/-  '  3_/  3U>. 

p\  3^  "  “  -  ^ou  v  1  + 


and,  for  the  five  second-order  quantities: 
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and  Re  is  the  turbulent  Reynolds  number  q  /v ,  with  being  related  to  the 
dissipation  rate  by 


3  2  3 

$  =  0.4  +  5itv  ^  =  0.4  (l  +  12.5  it/Re  )  (8) 

It  may  be  noted  that,  in  the  limit  Re  -►  «  ,  Eqs.  (1)  through  (7)  are 
somewhat  similar  to  those  derived  by  Launder  et  al.6  7  However,  low  Reynolds 
number  effects  are  also  included  here,  in  the  molecular  diffusion  terms  and  in 
the  Re  dependence  of  the  sink  terms  involving  CE  and  C<j> .  The  term 
14.8  pvq2/(<fy2)  is  a  "wall  term, "  required  to  obtain  a  well-behaved  solution 
in  the  viscous  sublayer,  y  -*•  0. 


For  the  high-speed  flows,  it  is  also  necessary  to  describe  the  tempera¬ 
ture  or  enthalpy  field.  The  mean  temperature  is  required  to  determine  density 
and  the  Reynolds  heat  flux  v ' T '  is  of  primary  interest  for  heat  transfer  con¬ 
siderations.  To  accomplish  this,  we  include  as  additional  dependent  variables 
the  mean  static  enthalpy  h,  the  mean-square  fluctuating  enthalpy  h ' 2 ,  and  the 
transverse  and  axial  components  of  the  Reynolds  heat  flux  v' h‘  and  u'h' . 
Required  closure  approximations  have  been  carried  out  in  a  manner  analogous  to 
those  leading  to  the  velocity  equations  above,  although  the  paucity  of  meas¬ 
urements  of  fluctuating  temperatures  has  made  it  difficult  to  completely 
verify  closure  approximations  against  basic  laboratory  data.  The  resulting 
thalpy  equations  are : 
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1.32  +  7.5  ir/Re.  1.165  +  12.5  ir/Re. 

_  _ A  _ A^ 

where  C  -  y  +  12  5  n/Re^  »  CT  "  1  +  12.5  Tt/ReA 

It  should  be  noted  that  terms  involving  fluctuating  densities  ( p ' )  have 
been  dropped  in  deriving  Eqs.  (1)  through  (12).  This  is  generally  permissible 
if  the  edge  Mach  number  is  below  4  or  5.  However,  if  need  be,  the  dominant 
effects  of  density  fluctuations  can  be  included  by  defining  a  generalized 
Reynolds  stress  R^j  =  pu^u  j/p  =  ujuj  +  p  '  uju  j/p .  Once  this  is  done,  the  pri¬ 
mary  effect  of  density  fluctuations  is  contained  in  a  relatively  unimportant 
diffusional  term  involving  p 1 v ' ,  which  can  be  related  to  v  T  .  The  resulting 
formulation  yielded  good  comparisons8  with  the  measurements  of  Hcrstman  et 
al.9  in  a  boundary  layer  at  Me  =  7. 

Boundary  conditions  to  Eqs.  (1)  through  (12)  are  generally  obvious: 


fluctuating  quantities  are  zero  at  a  sold  wall  or  at  the  outer  edge  (if  there 
is  no  free-stream  turbulence). 


3.  INVERSE  BOUNDARY  LAYER  CALCULATIONS 


The  inverse  boundary  approach  offers  an  attractive  approach  for  solving 
separated  boundary  layer  flows.  For  most  hypersonic  flows,  separation  zones 
are  relatively  thin.  Normal  pressure  variations  are  small,  and  boundary  layer 
theory  should  provide  a  useful  approximation.  Reyhner  and  Flugge-Lotz1 0 
pioneered  the  use  of  boundary  layer  methods  for  separated  regions,  by 
neglecting  streamwise  convection  terms  in  reverse  flow  zones. 

The  major  challenge  to  calculating  separated  flows  is  to  couple  the 
inner,  viscous  solution  to  the  outer,  inviscid  solution.  The  conventional 
direct  boundary  layer  approach  (pressure  or  edge  velocity  specified  as  a  func¬ 
tion  of  streamwise  distance)  is  singular  at  a  separation  point.  This  singu¬ 
larity  precludes  coupling  the  viscous  and  inviscid  solution.  However,  as 
originally  pointed  out  by  Catherall  and  and  Mangier,11  an  inverse  calculation 
(wherein  the  boundary  layer  thickness  or  wall  shear  is  prescribed  vs.  dis¬ 
tance)  is  regular  at  separation  points.  Thus,  it  should  be  possible  to  couple 
the  viscous  and  inviscid  solutions.  Of  course,  the  displacement  thickness  6* 
would  be  the  appropriate  boundary  layer  thickness  to  specify  for  this  coupling 
process . 

Wigton  and  Holt12  have  outlined  how  this  coupling  should  be  carried  out. 
One  first  computes  the  inverse  boundary  layer  flow  given  an  approximate  6* 
consistent  with  the  flow  configuration.  Then  one  computes  the  inviscid  flow 
past  the  shape  that  includes  the  same  5*.  In  general,  these  two  calculations 
will  yield  different  pressure  distributions.  Then,  based  on  some  knowledge  of 
the  viscous  and  inviscid  flow  characteristics,  a  revised  5*  is  selected  to 
yield  a  better  pressure  match.  The  iteration  continues  until  both  viscous  and 
inviscid  solutions  yield  essentially  the  same  pressure  fields. 

Several  inverse  boundary  layer  methods  have  been  developed  in  recent 
years.  Cebeci13  use<j  an  inverse  approach  to  compute  attached  compressible 
boundary  layers,  both  laminar  and  turbulent,  with  either  Cf(x)  or  5*(x) 


prescribed.  A  Newton's  method  was  used  to  determine  the  unknown  edge  pres¬ 
sure  at  each  station,  which  involved  inverting  a  matrix  of  derivatives  of  the 
finite-difference  equations  with  respect  to  the  unknown  pressure.  While  quite 
sound,  this  matrix  inversion  would  be  cumbersome  to  apply  to  the  present 
Reynolds  stress  formulation.  It  would  be  necessary  to  differentiate  all  of 
our  equations  for  turbulent  variables  with  respect  to  the  pressure;  this  would 
greatly  expand  the  required  matrix  inversion  process,  which  is  already  compli¬ 
cated  by  the  stiff  nature  of  the  governing  equations  for  the  Reynolds  stress 
variables . 

Cebeci,  Keller,  and  Williams1 ^  (see  also  Cebeci  and  Stewartson1 5 ) 
obtained  inverse  boundary  layer  solutions  of  separated,  incompressible  laminar 
flows.  The  velocity  components  are  transformed  to  the  standard  stream  func¬ 
tion  as  the  dependent  variable.  Specification  of  the  displacement  thickness 
6*(x)  yields  an  outer  boundary  condition  on  the  stream  function: 

'M  x ,  ye )  -  pUe(x)(ye  -  6*(x))  for  y  *  ye  .  (13) 

This  follows  directly  from  the  definitions  of  stream  function  and  displacement 
thickness,  and  provides  an  additional  boundary  condition  to  the  usual  ones 
(U»V=0  at  y=0,  U=Ue  at  y=»ye).  Thus  the  boundary  conditions  are  overspecified, 
and  the  edge  velocity  Ue(x)  is  obtained  as  an  eigenvalue.  Cebeci,  Keller  and 
Williams  showed  good  agreement  with  earlier  results  of  Carter  (to  be  discussed 
shortly),  with  far  fewer  iterations  than  required  by  Carter's  methods. 

While  quite  sound,  the  Cebeci,  Keller  and  Williams  method  is  not  easily 
extended  to  compressible  flows.  Equation  (13)  still  holds,  with  p  replaced  by 
pe(x).  But,  use  of  the  streamfunction  as  the  dependent  variable  results  in 
the  appearance  of  partial  derivatives  of  the  density  in  the  momentum  equation. 
These  derivatives  could  be  eliminated  by  use  of  a  density- transformed  y 
coordinate 


but  then  it  would  be  necessary  to  perform  an  integral  to  obtain  ye  in 
Eq.  (13): 


dY 

P„ 


Such  terms  involving  density  variations  are  quite  significant  in  hypersonic 
flows,  where  the  density  can  vary  by  large  factors  across  the  boundary  layer 
or  in  the  streamwise  distance.  Of  course,  an  energy  equation  must  be  solved 
simultaneously,  which  provides  the  required  thermodynamic  information  to 
determine  the  density.  However,  the  density  variations  make  it  difficult  to 
apply  the  inverse  calculations  with  this  method.  For  example,  suppose  the 
integral  of  Eq.  (15)  is  inserted  in  Eq.  (13)  for  ye.  The  variation  of  the 
integral  of  / p-1  dY  with  downstream  distance  can  be  as  large  as  changes  in  6*, 
at  least  over  certain  regions.  This  makes  it  very  difficult  in  practice  to 
obtain  reliable,  accurate  inverse  solutions  for  high-speed  flows  with  such  an 
approach . 


In  a  series  of  papers,  Carter  (see  Ref.  16)  has  pursued  a  somewhat 
different  approach  to  inverse  calculations,  which  he  has  successfully  extended 
to  compressible  laminar  and  turbulent  (with  an  eddy  viscosity  model)  applica¬ 
tions.  Carter  uses  the  density-transformed  normal  coordinate  (Eq.  (14))  and 
the  standard  streamfunction 

pU;  —  =  -PV  .  (16) 

With  these  definitions,  the  convective  operator  becomes 
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Carter  then  introduces  a  perturbation  streamfunction: 


Pe0[/  f-  dY  -  6  ] 
o 


mmsm 


or ,  rear rang ing , 


-  Y 


-  O 


(19) 


n  -  */«*.  t*  -  /  -  i)  dn  . 

o  p 


The  perturbation  streamf unction  obeys  the  following  differential  equation: 


a;  *, 

—  «-pU6(n-1  +  t)  — - - 
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(20) 


and  satisfies  the  outer  boundary  condition 


¥  0  as  n  +  00  •  (21) 

The  major  advantage  of  the  perturbation  streamf unction  is  that  the  pre¬ 
scribed  displacement  thickness  6*(x)  is  automatically  achieved;  this  may  be 
seen  by  comparing  Eqs.  (18)  and  (21)  with  Eqs.  (13)  and  (15)  above.  Thus  the 
specified  boundary  layer  growth  is  basically  satisfied  through  a  boundary  con¬ 
dition  (Eq.  (21)),  at  the  expense  of  an  additional  first-order  partial  differ¬ 
ential  equation  (Eq.  (20)).  Also,  Eq.  (17)  requires  3f/3x,  which  may  be 

A 

related  to  "P  and  6*  through  Eq.  (19). 


In  the  initial  stages  of  this  program,  we  implemented  Carter's  method  for 
laminar  boundary  layers  and  reproduced  some  of  his  earlier  results.  Figures  1 
and  2  show  calculations  for  a  separated  laminar  boundary  layer,  compared  to 
Carter's  solution  and  the  Navier-Stokes  computation  by  Briley.17  Not  only  was 
Carter's  calculation  reproduced,  but  also  the  good  agreement  with  the  more 
exact  Navier-Stokes  solution  indicates  the  utility  of  the  boundary  layer 
approach  for  flows  with  thin  separation  zones. 


As  already  mentioned,  Carter  extended  his  method  to  turbulent  flows, 
using  an  eddy  viscosity  formulation.  However,  we  found  the  method  to  be 
difficult  to  apply  to  the  more  sophisticated  Reynolds  stress  turbulence 
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Figure  2.  Velocity  profile  comparison  a 


formulation.  The  complications  arise  from  some  of  the  numerical  constraints 
associated  with  the  solution  of  the  five  simultaneous  turbulence  variables. 
These  equations  are  nonlinear,  strongly  coupled,  and  often  "stiff"  (i.e.,  the 
right  sides  are  dominated  by  production  and  dissipation  terms,  which  nearly 
cancel) . 


The  normal  velocity  convective  term  involves  S^/Sx,  which  is  related  to 
the  perturbation  streamfunction  by 
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In  solving  the  Reynolds  stress  equation,  it  is  frequently  necessary  to  take 
very  small  steps  in  the  streamwise  direction.  This  is  required  for  successful 
use  of  Newton's  method  for  solving  the  turbulence  equations.  Typically,  the 
step  size  Ax  is  chosen  so  that  none  of  the  dependent  variables  changes  by  more 
than  a  specified  percentage  (e.g. ,  10  or  20  percent)  from  one  step  to  the 
next.  In  the  start-up  process,  the  initial  step  sizes  are  no  larger  than  the 
smallest  Ay  values,  to  allow  a  smooth  adjustment  of  the  unknown  values  to 
their  "stiff"  solutions.  For  a  high  Reynolds  number  turbulent  boundary  layer, 
the  smallest  Ay  in  the  viscous  sublayer  can  be  several  orders  of  magnitude 
less  than  the  boundary  layer  thickness.  Of  course,  the  step  size  Ax  subse¬ 
quently  grows  exponentially  to  reasonably  large  values.  But,  when  taking 
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small  steps,  the  terms  3 *F/3 x  and  3U/3x  in  Eq.  (22)  can  have  quite  large 
values.  Physically,  3*F/3x  should  be  well  behaved,  which  means  that  the  right 
side  of  Eq.  (22)  can  be  expected  to  involve  small  differences  of  large  num¬ 
bers.  Thus,  with  small  Ax  steps  the  convective  operator,  the  left  side  of  the 
governing  equations,  is  stiff .  We  were  not  successful  in  developing  a  reli¬ 
able  numerical  method  for  solving  such  a  "doubly-stif f "  set  of  equations.  The 
difficulty  does  not  arise  if  small  steps  are  avoided,  as  Carter  was  able  to 
do  with  an  eddy  viscosity  formulation.  However,  the  perturbation  stream- 
function  approach  implies  a  minimum  acceptable  step  size,  whereas  the  Reynolds 
stress  equations  require  a  maximum  step  size  value,  and  these  two  requirements 
are  incompatible. 
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A  streamfunction  formulation  was  developed  to  avoid  the  difficult 
numerical  aspects  of  Carter's  perturbation  streamfunction,  as  well  as  the  com¬ 
plications  associated  with  use  of  the  Cebeci,  Keller  and  Williams1 4  boundary 
condition,  Eq.  (13).  A  normalized  incompressible  normal  coordinate  is  used: 
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and  the  normal  velocity  is  eliminated  in  favor  of  the  standard  streamfunction 
The  convective  operator  becomes 
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If  the  dependent  variables  and  density  are  normalized  by  pe(x)  and  Ue(x) 


U  -  U/U  (x) 


(25) 


p  =  p/p  (x) 
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etc. 


then  the  governing  equations  for  velocity  and  turbulent  kinetic  energy  take 
the  following  form 
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The  corresponding  equations  for  the  other  Reynolds  stress  variables  and  for 
the  enthalpy  unknowns  take  similar  forms.  In  addition,  since  the  stream- 
function  is  included  as  an  unknown  variable  in  Eqs.  (26)  and  (27),  an  addi¬ 
tional  partial  differential  equation  must  be  included,  which  is  simply  derived 
from  its  definition: 

p  U  5  -  U  (28) 

e  e  3n 

t  -  0  at  n  *  0  . 

This  formulation,  as  explained  thus  far,  does  not  provide  a  method  for 
performing  an  inverse  calculation.  Both  Ue(x)  and  5*(x)  are  included  as 
coefficients.  If  6*(x)  is  to  be  specified,  then  a  technique  must  be  provided 
for  determining  Ue(x).  To  this  end,  a  rather  straightforward  approach  was 
taken,  employing  what  is  basically  a  secant  method.  At  each  new  station,  the 
solutions  are  obtained  at  successive  values  of  Ue  in  an  iterative  fashion 
until  the  prescribed  6*  is  obtained.  With  a  secant  method,  two  values  of  ue 
are  selected  at  the  new  (unknown)  station.  The  governing  equations  are  solved 
for  each  value,  and  the  resulting  5*  values  are  computed.  From  these  two 
values,  one  interpolates  to  a  new  Ue  value  corresponding  to  the  desired  6*- 
The  equations  are  solved  again  for  this  new  5*,  and  subsequent  iterations  per¬ 
formed  until  the  prescribed  5*  value  is  achieved  to  an  acceptable  tolerance. 

In  practice,  after  the  first  step,  it  is  possible  to  retain  the  slope  of  the 
curve  in  Figure  3  from  one  step  to  the  next,  so  as  to  avoid  wasting  two 
initial  solutions  before  interpolating  to  the  desired  5*.  Note  that  for  very 
small  steps  when  the  turbulence  variables  are  changing  rapidly,  5*  and  Ue 
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Figure  3.  Sketch  of  secant  method  for  inverse  boundary  layer  solution 

change  only  slightly  from  one  step  to  the  next  and  the  secant  method  converges 
quickly.  Note  also  that  this  approach  is  a  straightforward  (and  perhaps  less 
elegant)  alternative  to  Cebeci's  method.13  In  essence,  he  obtained  the  slope 
on  Figure  3  by  an  elaborate  matrix  inversion.  His  technique  might  be  more 
efficient  numerically  than  that  used  there,  but  would  certainly  entail 
significantly  greater  analytical  and  programming  effort. 

In  Figure  4  we  compare  the  direct  and  inverse  solutions  for  a  Mach  7  cold 
wall  boundary  layer.  The  direct  calculation  was  first  performed  for  a  flat 
plate,  with  no  pressure  gradient.  Using  the  displacement  thickness  from  this 
direct  calculation  as  input,  the  inverse  solution  was  then  obtained  according  i 

to  the  method  just  outlined.  The  skin  friction  coefficients  are  in  good 
agreement,  but  the  momentum  thickness  values  differed  by  10  to  20  percent.  Of 
course,  the  displacement  thickness  is  forced  to  be  identical  for  the  two 

i 

calculations,  so  the  shape  factors  are  different.  This  difference  undoubtedly 
reflects  the  numerical  accuracy  of  these  calculations,  which  involved  rela¬ 
tively  coarse  grids  (39  points  across  the  boundary  layer).  A  finer  mesh  1 

appears  to  be  necessary  for  accurate  values  of  the  momentum  thickness.  j 
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4.  ADVERSE  PRESSURE  GRADIENTS  EFFECTS 


A  major  issue  that  seems  particularly  critical  is  the  ability  of  present 
models  to  describe  the  detailed  development  of  boundary  layer  turbulence  in 
the  presence  of  strong  adverse  pressure  gradients.  Pressure  gradients  cause 
significant  redistributions  of  energy  between  the  various  components  of  the 
fluctuating  turbulent  energy.  The  "five-equation"  Reynolds  stress  model  used 
in  this  study  should  be  superior  to  simpler  models  and,  hopefully,  should  pro¬ 
vide  an  accurate  description  of  turbulence  behavior  in  separating  boundary 
layers.  The  closure  model  obeys  required  tensor  rules,  with  the  closure  con¬ 
stants  based  on  data  from  experiments  on  grid  turbulence  in  the  presence  of 
mean  shear  and  strain.  Figure  5  is  an  example  of  a  comparison  between  the 
model  and  measurements  on  grid  turbulence  subjected  to  plane  strain.1®  We 
also  examined  experimental  results  on  low-speed  boundary  layers  in  favorable 
pressure  gradients,  including  the  "relaminarization"  effect,  with  reasonable 
success.  However,  in  view  of  the  importance  of  the  turbulence  behavior,  we 
have  examined  several  laboratory  experiments  on  attached  boundary  layers  in 
adverse  pressure  gradients . 

Bradshaw  Case  C 


Bradshaw's  experiment19  consisted  of  a  low-speed  boundary  layer  on  a  flat 
plate,  artifically  tripped  well  upstream  of  the  region  of  interest.  Reynolds 
numbers  were  in  the  range  Re  -  10,000  to  20,000,  so  the  turbulent  boundary 
layer  should  be  well  developed.  The  inviscid  velocity  decreased  from 
128  ft/sec  to  about  90  ft/sec  over  a  distance  of  60  in.  Figure  6  shows  the 
variation  of  skin  friction  and  displacement  thickness  over  the  range  of  the 
adverse  pressure  gradient,  and  Figures  7  through  9  show  the  profiles  of  mean 
velocity,  Reynolds  stress  and  fluctuating  velocity  components  at  the  84  in. 
station,  at  the  end  of  the  adverse  pressure  gradient.  The  most  dissappointing 
result  is  the  computed  Reynolds  stress,  which  seems  to  be  "amplified"  too  much 
by  the  pressure  gradient.  Correspondingly,  the  computed  wall  shear  is  high, 
the  boundary  layer  becomes  thicker  than  observed,  and  the  mean  velocity 
profile  shape  is  not  in  detailed  agreement  with  Bradshaw's  measurements. 
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Figure  5.  Decay  of  grid  turbulence  subjected  to  plane  strain,  comparisons 
with  the  experiment  of  Tucker  and  Reynolds18 
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Samuel  and  Joubert  Experiment 

Another  experiment  on  an  attached  incompressible  boundary  layer  in  an 
adverse  pressure  was  conducted  by  Samuel  and  Joubert.20  In  Figure  10  we  show 
the  variation  of  skin  friction  and  displacement  thickness.  Profiles  of  mean 
velocity  and  Reynolds  stress  at  three  stations  (at  the  onset,  in  the  middle, 
and  near  the  end  of  the  adverse  pressure  gradient)  are  shown  in  Figures  11 
and  12.  Here  again,  the  turbulent  shear  is  significantly  over-predicted  near 
the  end  of  the  adverse  pressure  gradient.  However,  the  computed  boundary 
layer  thickness  and  velocity  profiles  are  in  better  agreement  with  the 


measured  values. 

Lewis,  Gran,  and  Kubota 

This  experiment21  involved  an  attached  turbulent  boundary  layer  in  an 
adverse  pressure  gradient,  at  supersonic  speeds.  The  boundary  layer  edge  Mach 
number  decreased  from  4.0  to  3.0,  followed  by  a  region  of  favorable  pressure 
gradient.  Figures  13  through  15  show  the  streamwise  variations  of  skin 
friction,  displacement  thickness  and  momentum  thickness.  In  this  case,  in  the 
more  downstream  region,  the  skin  friction  is  underpredicted  and  the  boundary 
layer  thicknesses  are  overpredicted. 
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Figure  12.  Reynolds  stress  profiles  for  Samuel  and  Joubert  experiment 
at  X  *  3.4  m 


5.  ATTEMPTS  AT  CALCULATING  SEPARATED  FLOWS 


The  most  difficult  aspect  of  this  work,  as  one  might  expect,  is  the 
application  of  the  inverse  method  to  separated  turbulent  flows.  As  a  first 
test  case,  we  selected  an  incompressible,  adiabatic  case  at  a  Reynolds  number 
of  106,  as  calculated  already  by  Carter.16  Although  no  experimental  data  are 
available  at  the  precise  conditions,  reproducing  Carter's  results  with  the 
rather  different  present  model  would  seem  to  be  a  useful  goal. 

Carter's  published  displacement  thickness  was  specified  as  input 
(Figure  16).  The  initial  boundary  layer  was  taken  to  be  a  classical  flat 
plate  boundary  layer  developed  in  the  absence  of  a  pressure  gradient.  The 
surprising  result  of  our  inverse  calculations  was  that  the  present  model  does 
not  predict  separation  at  any  point  during  the  substantial  increase  of 
boundary  layer  thickness.  Figure  17  shows  the  computed  skin  friction,  com¬ 
pared  to  Carter's  published  values.  Whereas  the  solution  obtained  by  Carter 
showed  separation  at  an  X  value  slightly  greater  than  1.1,  no  such  result 
occurs  with  our  calculation.  In  Figure  18  we  show  the  edge  velocity  derived 
from  the  solution.  For  the  present  solution,  a  rather  more  gradual  adverse 
pressure  gradient  is  implied. 

The  significant  characteristic  of  our  solution  in  the  region  where 
separation  would  be  expected  is  the  behavior  of  the  mean  velocity  profile  and 
of  the  Reynolds  stress.  Figure  19  shows  the  mean  velocity  profiles  at  two 
selected  stations,  along  with  the  corresponding  Reynolds  stresses  in 
Figure  20.  The  relative  Reynolds  stress  is  amplified  by  the  adverse  pressure 
gradients.  At  the  same  time,  the  mean  velocity  develops  steep  gradients  at 
the  wall  (supporting  the  sustained  high  shear)  and  a  "plateau"  with  even  a 
mild  "valley". 

The  mean  velocity  and  Reynolds  stress  development  are  certainly  coupled 
in  a  complex  manner,  and  could  easily  be  an  anomolous  consequence  of  the 
Reynolds  stress  model  employed  here.  The  development  of  the  Reynolds  stress 
is  reminescent  of  the  behavior  in  the  attached  adverse  pressure  gradient 
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experiments  of  Bradshaw19  and  Samuel  and  Joubert,20  discussed  above.  An 
obvious  issue  is  the  sensitivity  of  these  results  to  the  turbulence  modeling 
closure  assumptions  and  approximations. 


The  most  surprising  aspect  of  our  turbulence  solutions  in  adverse  pres¬ 
sure  gradients  is  the  amplification  of  the  Reynolds  stress  that  was  discussed 
above  in  connection  with  the  Bradshaw  and  Samuel  and  Joubert  experiments,  and 
with  the  present  attempt  to  reproduce  Carter's  separated  flow  solution.  This 
suggests  a  shortcoming  in  the  turbulence  closure  modeling.  However,  it  is  not 
obvious  where  the  difficulty  lies.  Note  that  there  is  no  direct  pressure 
gradient  OU/3x)  term  in  the  equation  for  the  Reynolds  stress  equation  (u  v  ) , 
in  contrast  to  the  equations  for  the  components  of  the  turbulent  kinetic 
energy  and  the  dissipation  function.  This  is  a  direct  result  of  the  tensor 
form  assumed  for  the  so-called  pressure-strain  terms  (same  form  introduced  by 
Hanjalic  and  Launder6).  One  could  arbitrarily  introduce  a  term  in  the  equation 
governing  u*  v'  to  damp  turbulence  in  regions  of  adverse  pressure  gradients 
(i.e.,  a  term  that  would  reduce  the  magnitude  of  u  v  when  3U/3x  is  nega¬ 
tive).  But  that  would  violate  the  desired  tensor  character  of  the  pressure- 
strain  term  and  there  is  no  rational  basis  for  such  as  approach.  The  primary 
role  of  such  terms  is  dictated  by  kinematic  considerations,  and  can  be  derived 
without  recourse  to  "closure"  approximations.  We  investigated  the  effect  of 
variations  in  closure  constants  appearing  in  other  terms  in  the  five-equation 
model.  Changing  the  wall  term  in  the  dissipation  equation,  or  the  turbulent 
diffusion  terms,  does  not  resolve  the  discrepancies,  and  would  jeopardize 
agreement  with  other,  more  basic  data  upon  which  the  closure  constants  are 
based.  Arbitrarily  increasing  the  3U/3x  term  in  the  dissipation  equation 
removes  the  anomolous  velocity  profile  shown  in  Figure  19,  but  this  again 
violates  a  required  tensor  variation.  However,  none  of  our  attempts  at  vary¬ 
ing  closure  approximations  produced  results  that  are  qualitatively  similar  to 
those  obtained  by  Carter, 16  or  even  produced  a  well-defined  separation  zone. 


6.  RECOMMENDATIONS  FOR  FUTURE  WORK 


It  should  be  clear  from  the  results  of  this  study  that  the  implementation 
of  an  inverse  boundary  layer  method  for  high-speed  turbulent  separated  flows 
is  a  more  ambitious  goal  than  could  be  accomplished  in  this  program.  We 
firmly  believe  that  an  inverse  method  can  be  developed  as  a  useful  engineering 
tool.  However,  three  key  areas  remain  where  advances  are  clearly  required. 

Turbulence  Model  -  The  type  of  Reynolds  stress  model  employed  here  should 
provide  a  good  description  of  turbulence  behavior  in  the  presence  of  pressure 
gradients.  The  use  of  a  five-equation  model  should  account  for  the  influence 
of  pressure  gradients  on  the  degree  of  isotropy,  the  distribution  of  energy 
between  the  different  velocity  components,  and  the  Reynolds  stress;  the 
simpler  two-equation  (k-e)  approaches  inherently  assume  an  unvarying  degree  of 
isotropy  and  should  be  less  accurate.  However,  comparisons  of  results  of  the 
present  five-equation  formulation  with  measurements  on  attached  boundary 
layers  in  adverse  pressure  gradients  are  very  dissappointing.  More  experi¬ 
ments,  particularly  on  supersonic  boundary  layers,  are  clearly  in  order. 

Also,  some  rethinking  of  the  basic  formulation  of  current  Reynolds  stress 
modeling  is  needed.  For  example,  it  may  be  necessary  to  account  for 
anisotropic  dissipation  effects  and/or  varying  turbulence  length  scales  in 
different  directions  in  the  presence  of  strong  pressure  gradients. 

Inverse  Boundary  Layer  Methods  -  A  truly  acceptable  method  for  performing 
inverse  boundary  layer  computations  has  yet  to  be  developed.  The  current 
methods  (at  least  for  specified  displacement  thickess)  require  several  times 
as  much  computational  effort  as  the  corresponding  direct  calculation.  More 
iterations  are  required,  and  there  are  various  step  size  constraints.  In  our 
experience,  the  inverse  computation  is  not  very  "robust". 

The  primary  complication  with  the  inverse  approach  lies  with  the 
displacement  thickness.  This  quantity  controls  the  interaction  between  the 
inner  viscous  region  and  the  outer  inviscid  region,  so  the  displacement 


thickness  is  the  input  of  interest  for  the  inverse  solution,  rather  than  skin 
friction  or  mass  flux  within  the  boundary  layer.  Unfortunately,  6*  is  defined 
by  an  integral  across  the  boundary  layer.  It  does  not  naturally  appear  as  a 
boundary  condition  or  as  a  coefficient  in  the  boundary  layer  equations. 
Carter's  approach  introduces  6*  into  the  pde's,  but  in  the  form  of  a  small 
difference  of  large  numbers.  Cebeci's12  formulation  involves  a  complex 
iterative  Newton's  method  (matrix  inversion),  and  our  method  requires  an  addi¬ 
tional  (brute  force)  iteration.  The  approach  of  Cebeci,  Keller  and  Williams14 
is  difficult  to  extend  to  compressible  flows.  Clearly,  an  improved  approach 
is  necessary  if  the  inverse  method  is  to  find  wide  application. 

Iteration  Between  Viscous/Inviscid  Flows  -  A  final  requirement  is  a 
reliable  algorithm  for  matching  the  boundary  layer  solution  to  the  outer 
inviscid  flow.  Existing  coupled  solutions  to  date  have  been  for  transonic 
flows,  and  have  been  obtained  by  trial  and  error  iterations.  Wigton  and 
Holt1 2  have  established  the  theoretical  foundation  for  a  method  to  couple  the 
inner  and  outer  flows,  but  no  algorithms  have  been  developed  for  hypersonic 
flows.  We  had  so  much  difficulty  with  the  inverse  method  that  it  was  not 
possible  to  even  address  this  issue  in  this  effort.  Since  the  Newtonian 
approximation  (pressure  proportional  to  the  square  of  the  sine  of  the  local 
body  angle)  should  be  accurate  for  most  hypersonic  conditions,  it  is  tempting 
to  expect  that  an  algorithm  based  on  the  local  derivative  of  6*  can  be  devel¬ 
oped.  On  the  other  hand,  upstream  influence  transmitted  through  the  separa¬ 
tion  zone  may  necessitate  a  more  complex  algorithm.  If  a  relatively  straight¬ 
forward  algorithm  is  not  possible,  then  boundary  layer  methods  will  not  prove 
to  be  useful  for  thin  separated  flows. 

Our  experience  with  the  coupling  schemes  found  in  the  literature 
(Carter16,  Wigton  and  Holt1 2)  was  not  very  satisfying  for  our  hypersonic 
boundary- layer  interactions.  The  initial  attempts  at  utilizing  the  rational 
Wigton  and  Holt  scheme  have  proved  to  be  inadequate.  Wigton22  himself 
abandoned  his  scheme  in  favor  of  Carter's^  simpler  method  when  solving  three- 
dimensional  transonic-wing  viscous-inviscid  interaction  problems.  These 
calculations  necessitated  hundreds  of  iterations  to  convergence;  some  even 


took  103  cycles  to  obtain  acceptable  answers.  The  methods  indicated  above  are 
local  coupling  schemes.  What  is  needed  to  tame  the  iteration  cycles  is  a 
global  scheme  wherein  changes  in  the  displacement  thickness  (or  boundary-layer 
mass  transpiration  puS*)  influence  changes  along  the  entire  boundary  layer. 
Local  coupling  schemes  may  encounter  difficulties  due  to  their  inability  to 
enforce  a  downstream  boundary  condition  (Wigton23).  Werle  and  Vatsa24  have 
anticipated  these  difficulties  and  specified  a  downstream  boundary  condition 
in  order  to  compute  interacting  supersonic  laminar  boundary-layer  separations 
with  a  time -dependent  ADI  technique.  New  coupling  algorithms  must  face  the 
global  nature  of  these  interaction  problems  before  inverse  methods  will  prove 
to  be  useful  for  separated  turbulent  flows. 
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